Spherical coupled-cluster theory for open-shell nuclei 
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Background A microscopic description of nuclei is important to understand the nuclear shell-model from fundamental princi- 
ples. This cannot be achieved without an effective approximation scheme, especially if three-nucleon forces are included. 

Purpose Define and evaluate an approximation scheme that can be used to study nuclei that are described as two particles 
attached to a closed (sub-)shell nucleus. 

Methods The equation-of-motion coupled-cluster formalism has been used to obtain ground and excited state energies. This 
method is based on the diagonalization of a non-Hermitian matrix obtained from a similarity transformation of the 
many-body nuclear Hamiltonian. A chiral interaction at the next-to-next-to-next-to leading order (N'^LO) using a cutoff 
at 500 MeV was used. 

Results The ground state energies of ®Li and ®He were in good agreement with a no-core shell-model calculation using the 
same interaction. Several excited states were also produced with overall good agreement. Only the — 3^ excited state 
in Li showed a sizeable deviation. The ground state energies of 1*0, and i«Ne were converged, but underbound 
compared to experiment. Also the calculated spectra were converged and comparable to both experiment and shell- model 
studies in this region. Some excited states in were high or missing in the spectrum. 

Conclusions Converged results were obtained for all nuclei and the method can be used to describe nuclear states with simple 
structure. Especially the ground state energies were very close to an exact diagonalization. To obtain a closer match 
with experimental data, effects of three-nucleon forces, the scattering continuum as well as additional configurations in 
the coupled-cluster approximations, are necessary. 

PACS numbers: 21.60.De,21.10.Dr,21.60.Gx,31.15.bw 
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I. INTRODUCTION 

In the past decade, the amount of computing resources 
made available for scientific research, has grown several 
orders of magnitude. This trend will continue during 
this decade, culminating in exascale computing facihties. 
This will promote new insights in every disciphne as new 
problems can be solved and old problems can be solved 
faster and to a higher precision. 

In nuclear physics, the long standing goal of a pre- 
dictive theory, seems to approach fulfillment ever faster. 
However, complex obstacles have yet to be overcome. An 
exact solution to the nuclear many-body problem, is only 
possible for light nuclei (See Leidemann and Orlandini [l] 
for a recent review on many-body methods). With a fi- 
nite basis expansion, a full diagonalization can currently 
be performed for nuclei in the p-shell region Q , with the 
possibility of approaching light sd-shell nuclei with suffi- 
cient computing resources. For ab initio access to larger 
nuclei, the problem has to be approached differently. 

In coupled-cluster theory, a series of controlled approx- 
imations are performed to generate a similarity transfor- 
mation of the nuclear Hamiltonian, where the ground 
state energy of a closed (sub-)shell reference nucleus is 
singled out. The similarity-transformed Hamiltonian is 
then diagonalized to calculate excited states and states 
of nuclei with one or more valence nucleons. This defines 
the equation-of-motion coupled-cluster (called EOM-CC 
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from now on) framework (see Bartlett and Musial [3] for 
a recent review and Shavitt and Bartlett [4] for a text- 
book presentation). Recently, this method was applied 
to the oxygen Q and the calcium [6] isotopic chains, ex- 
tending the reach of ab initio methods in the medium 
mass region. Calculations in the nickel and tin regions 
are also planned. 



In this work, we refine the EOM-CC method for two 
valence nucleons attached to a closed (sub-)shell refer- 
ence. The general theory was presented in Jansen et al. 
0, where calculations were performed in small model 
spaces. The working equations are now completely re- 
worked in the spherical formalism, where analytical opti- 
mizations due to the rotational invariance of the Hamil- 
tonian, enable calculations in large model spaces. This 
involves a fair bit of angular momentum algebra, where 
the potential for errors is quite large. We therefore explic- 
itly include all equations and transformations that have 
been used in the spherical formalism for future reference. 

Section |lT] gives a brief overview of general coupled- 
cluster theory and the equation-of-motion extensions, 
emphasizing the parts that are necessary to develop the 
spherical scheme presented in Section IIIII In Section IIVI 
we discuss numerical results obtained for selected states 
in A — 6 and A — 18 systems, emphasizing the conver- 
gence properties. Finally, in Section |V] we present our 
conclusions and discuss the road ahead. We also include 
an appendix, where angular momentum transformations 
used in Section Hill are defined. 
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II. COUPLED-CLUSTER THEORY 

In this section we introduce the Hamiltonian that en- 
ters our calculation, a brief review of single reference cou- 
pled cluster-theory together with the equation-of-motion 
(EOM-CC) extensions to coupled cluster theory that al- 
low the calculation of excited states. In this framework, 
also nuclei with different particle numbers can be ap- 
proached by choosing the appropriate Slater determinant 
basis. We keep the presentation short and concentrate on 
the aspects important for deriving the spherical version 
of the 2PA-E0M-CCSD method presented in section Iml 
The interested reader is referred to Shavitt and Bartlctt 
Q for details on coupled-cluster theory. 

We use the intrinsic Hamiltonian 



Throughout this work we will use the convention that 
the indices ijk . . . denote states below the Fermi level 
(holes) , while the indices abc . . . denote states above the 
Fermi level (particles). For an unspecified state, the in- 
dices pgr... are used. The amplitudes t1^"'ij^ will be 
determined by solving the coupled-cluster equations. In 
the singles and doubles approximation we truncate the 
cluster operator as 



T: 



CCSD 



Ti 



(6) 



which defines the coupled-cluster approach with singles 
and doubles excitations, the so-called CCSD approxima- 
tion. The unknown amplitudes result from the solution 
of the non-linear CCSD equations given by 



H 



A* J ^ 2m 



A 



mA* 



(1) 



0, 
. 



Here A is the number of nucleons in the reference state. 



The term 



A* is the mass number of the target nucleus, and v^ is jj = exp (— T)Hexp (T) = ( Hexp (T) 

flip nnrlpnTi-mir'lpnn infprpirfinn AA^p ^srill limif nnrQpKrPQ ^ 



the nucleon-nucleon interaction. We will limit ourselves 
to two-body interactions only. In second quantization, 
the Hamiltonian can be written as 



(7) 



(8) 



H = 



pq 



(2) 



The term (pgljrs) is a shorthand for the matrix elements 
(integrals) of the two-body part of the Hamiltonian of 
Eq- ©, 97'' and s represent various single-particle 
states while stands for the matrix elements of the one- 



body operator in Eq. 
erators like aj, and Or 



([T]). Finally, second quantized op- 
create and annihilate a nucleon in 



the state q and p, respectively. These operators fulfill the 
canonical anti-commutation relations. 



A. Single-reference coupled-cluster theory 

In single-reference coupled-cluster theory, the many- 
body ground-state j^E'o) is given by the exponential 
ansatz, 



is called the similarity-transform of the normal-ordered 
Hamiltonian, but ([8]) only defines a proper similarity 
transformation when the cluster operator Q is used 
without truncation. The state ' ) is a Slater deter- 
minant that differs from the reference |$o) by holes in 
the orbitals ij . . . and by particles in the orbitals ah ... . 
The subscript C indicates that only connected diagrams 
enter. 

Once the and tf^ amplitudes have been determined 
from Eq. (O , the correlated ground-state energy is given 
by 



(9) 



l^-o) =exp(T)|$o 



(3) 



Here, Eq denotes the vacuum expectation value with re- 
spect to the reference state. 

The CCSD approximation is a very inexpensive 
method to obtain the ground state energy of a nucleus. In 
most cases however, the accuracy is not satisfactory Q. 
The obvious solution would be to include triples exci- 
tations in Eq. ([B]) to define the CCSDT approximation. 
This leads to an additional set of non-linear equations 



Here, |<I>o) is the reference Slater determinant, where all 
states below the Fermi level are occupied and T is the 
cluster operator that generates correlations. The opera- 
tor T is expanded as a linear combination of particle-hole 
excitation operators 



0, 



(10) 



T = Ti + T2 



(4) 



where T„ is the n-particle-n-hole(nP-nH) excitation op- 
erator 



(5) 



that has to be solved consistently. Unfortunately, such 
a calculation is very computationally expensive [9] . The 
computational cost of CCSDT scales as ?^onf^, where Uo 
is the number of single-particle states occupied in the 
reference determinant and n^, are the number of unoccu- 
pied states. For comparison, the computational cost of 
the CCSD approximation scales as n^n^. 

Instead of solving the coupled-cluster equations ([7]) 
including triples excitations, we calculate a correction 
to the correlated ground state energy ([5]), using the A- 
CCSD(T) approach [l^, Here, the left-eigenvalue 
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problem using the CCSD similarity-transformed Hamil- 
tonian is solved, yielding a non-iterative triples correction 
to the ground state energy. The left- eigenvalue problem 
is given by 



($o|AH = S($o|A, 
where A is a de-excitation operator, 
A = i + Ai -f A2, 

and 

Ai =^A>^a|, 

ia 

A2 =^A^>6aaala]. 

ijah 



(11) 

(12) 

(13) 
(14) 



The unknown amplitudes \\ and A^^j, are the compo- 
nents of the left-eigenvector with the lowest eigenvalue 
in Eq. (fTTj) . Once found, the energy correction is given 
by 



ahc\ 
ijkl 



ijkabc 



^ei^iJk\(yN%)j^0). (15) 



ijk 



Here, F^p is the part of the normal-ordered one-body 
Hamiltonian that annihilates particles and creates holes. 
The energy denominator is defined as 



^ijk — fa fjj ^" fkk faa fbb fc' 



where fpp are the diagonal elements of the normal or- 
dered one-body Hamiltonian F. 

Using this approach, the ground state wave func- 
tion ([2]) and the similarity transformed Hamiltonian (|5]) 
are calculated using the CCSD approximation, while the 
ground state energy is given by 



Eacc = E( 



cc 



(17) 



This approximation has proven to give very accurate re- 
sults for closed (sub-) shell nuclei |12| . 



B. Equation-of-motion coupled-cluster(EOM-CC) 

In nuclear physics, the single reference coupled clus- 
ter method defined by the coupled-cluster equations©, 
is normally used to obtain the ground state energy of 
a closed (sub-)shell nucleus. While it is possible to 
apply the CC method to any reference determinant to 
obtain the energy of different states, the equation-of- 
motion coupled-cluster (EOM-CC) framework is usually 
employed for such endeavors. 

When the cluster operator (U) is used without trun- 
cation, Eq. ([S]) defines a similarity transformation. This 



guarantees that the eigenvalues of H are equivalent to the 
eigenvalues of the intrinsic Hamiltonian ([IJ and that the 
eigenvectors are connected by the transformation defined 
by Eq. ([3]). When the cluster operator (U) is truncated, 
as in the CCSD approximation, Eq. (|8]) is an approxima- 
tion to a similarity transformation. This means that the 
eigenvalues of H are now only approximations to the ex- 
act eigenvalues. This is the foundation of the EOM-CC 
approach. 

To simplify the equations and for effective calculations, 
we define the eigenvalue problem in the EOM-CC ap- 
proach to be that of the difference between a target state 
and the coupled cluster ground state Formally, a 

general state of the A-body nucleus is written 



= n 



(18) 



Here is an excitation operator that creates the 
state l^*^), when applied to the coupled cluster ground 
state l^'o). The label fi identifies the quantum num- 
bers(energy, angular momentum, . . . ) of the target state. 
The Schrodinger equations for the target state and the 
coupled-cluster ground state are written 



He*|$o) =^cce*|$o). 



(19) 
(20) 



Here i?^ is the energy of the target state and Ecc is the 
coupled-cluster reference energy in Eq. 

We can now multiply Eq. (HH) by e'^ and Eq. (gHl) by 
Ctf^e~'^ from the left and take the difference between the 



(16) 

two equations 



H,r!„ 



1*0) 



(21) 



where = Ef^—Ecc ^^nd we have used that f2^, T = 0. 
Finally, none of the unconnected terms in the evaluation 
of the commutator survive, resulting in 



(22) 



This operator equation can be posed as a matrix eigen- 
value problem where uj^ are the eigenvalues and the ma- 
trix elements of Cl^ are the components of the eigenvec- 
tors. The subscript C implies that only terms where fl 
and are connected by at least one contraction survive. 
In diagrammatic terms, this means that only connected 

diagrams appear in the operator product ^Hf2^^ . 

The similarity-transformed Hamiltonian ([8]) is a non- 
hermitian operator and is diagonalized by the scheme 
formulated by the Arnoldi al gor ithm ( for details see for 
example Golub and Van Loan [ISj ) . This scheme relies on 
the repeated application of the connected matrix vector 
product defined by Eq. ([2^ . A left-eigenvalue problem 
is solved to obtain the conjugate eigenvectors [ij, but 
this is beyond the scope of this article. 
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To find the explicit expressions for the connected ma- 
trix vector product, the excitation operator must be 
properly defined. When used for excited states of an 
A-body nucleus, the excitation operator in Eq. is 
parametrized in terms of nP-nH operators. We define 



17„ = R = 1 + Ri + R2 



(23) 



where 



= E C:C<---4„«.„---«n- (24) 



Zl ,...ln 

ai ,...a„ 



The unknown amplitudes r(sub and superscripts 
dropped) are the matrix elements of R, 



define the particle attached equation-of-motion coupled- 
cluster ^16)] (PA-EOM-CC) and the particle removed 
equation-of-motion coupled-cluster pTI] (PR-EOM-CC) 
approaches. These methods have been used successfully 
in quantum chemistry for some time (see Bartlett and 
Musial for a review), but have also recently been im- 
plemented for use in nuclear structure calculations fl^. 

In Jansen et al. we defined the 2PA-E0M-CCSD 
and 2PR-E0M-CCSD methods for two particle attached 
to and removed from a closed (sub-)shell nucleus. For 
this problem, the excitation operators were given by 



where 



R 



A±2 



R 



A±2 



.R 



A±2 



(31) 



ai...a„ _ 



(25) 



and can be grouped into a vector that solves the eigen- 
value problem in Eq. ([22]) . The matrix vector product is 
found by looking at each individual element 



(26) 



and find the appropriate expressions using a diagram- 
matic approach. 

Calculations using the full excitation operator (|23|) , are 
not computationally tractable, so an additional level of 
approximation is introduced by a truncation. When the 
CCSD approximation is used to obtain the reference wave 
function, the excitation operator is truncated at the 2P- 
2H level [Tg and defines EOM-CCSD. 

In the EOM-CC approach, the states of A ± fc nuclei 
are also treated as excited states of an A-body nucleus. 
We write 



I* 



A±fc- 



(27) 



The operator fi^ and the energies E'^ of the target state, 
also solve the eigenvalue problem in Eq. ((22|) . The energy 



difference = i?^ — Eq is now the excitation energy 
of the target state in the nucleus A ± k, with respects 
to the closed-shell reference nucleus with the mass shift 
A* = A ± k in the Hamiltonian ([ij. This mass shift 
ensures the correct kinetic energy of the centcr-of-mass 
is used in computing the A± k nuclei. 
The operators 



R 



A±l 



. .R 



A±l 



(28) 



where 



R^+i = 



1 V r"!-"" 



(n!)(n - 1) 



ji,...i„_i 

ai , . ..an 



4. 



R: 



A-l 



{nl){n - 1) 



I \ " Qi...a„_i 

2-1)1 



ll ,...tn 
CLl ,.. .art — l 



(29) 



(30) 



R: 



A+2 



1 



(n!)(n-2) 



E ai...a„ 
!l...l„_; 



'il,...I„-2 



X 4i •■•aLa»„ 



R; 



A-2 



{n\){n- 2)1 



E 



ai...a„_2 
■ Zl . ..in 



r, 



2l ,. . .'In 

ai,...an-2 



X ai . . . al aj 



(32) 



(33) 



In this article, our focus is on the 2PA-E0M-CCSD 
method, where we truncate the operator at the 3P-1H 
level. This approximation is suitable for states with a 
dominant 2P structure. It is already computationally 
intensive with up to 10^ basis states for the largest nu- 
clei attempted. A full inclusion of 4P-2H amplitudes is 
therefore not feasible at this time. 



C. Spherical coupled cluster theory 

For nuclei with closed (sub-)shell structure, the refer- 
ence state has good spherical symmetry and zero total 
angular momentum. For these systems, the cluster oper- 
ator Q is a scalar under rotation and depends only on 
reduced amplitudes. Thus, 

Tl=E^^('^) WaiJ)®h{J)t 



and 



T, 



E 

ijabj 



ttiJ) 



4(ja)® 4(ifc) 



(g) [a,j{jj)(g)a,i{ji)Y 



where the amplitudes i(J)(sub- and superscript dropped) 
are a short form of the reduced matrix elements of the 
cluster operator (j4|) (See appendix lAl for details). Also, 
J is a label specifying the total angular momentum of a 
many-body state and standard tensor notation has been 
used to specify the tensor couplings. 

As the similarity-transformed Hamiltonian ([5]) is a 
product of three scalar operators (remember that the ex- 
ponential of an operator is defined in terms of its Taylor 
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expansion), it is also a scalar under rotation. This al- 
lows a formulation of the coupled-cluster equations that 
is completely devoid of magnetic quantum numbers, thus 
reducing the size of the single-particle space and the num- 
ber of coupled non-linear equations to solve in Eq. ([7]). 
For further details, see Hagen et al. jl2| . 

Within the same formalism, we wish to establish the 
connected operator product in Eq. ([^ . This will greatly 
reduce the computational cost of calculating the prod- 
uct, but also allow a major reduction in both the single- 
particle basis and the number of allowed configurations 
in the many-body basis. 

Given a target state with total angular momentum J 
(in units of he), the excitation operator, 51^ ([IH]), is a 
spherical tensor operator by definition(see for example 
Bohr and Mottelson [ISJ)- It has a rank of J, with 2J+1 
components labelled by the magnetic quantum number 
Me [-J, J]. We write 

n, = R^^>^{J,M), (34) 

where A is the number of particles in the reference state, 
AzLk is the number of particles in the target state, while fi 
identifies a specific set of quantum numbers. Identifying 
the excitation operator as a spherical tensor operator, 
invokes an extensive machinery of angular momentum 
algebra with important theorems. Of special importance 
is the Wigner-Eckart theorem(see for example Edmunds 
[20|). which states that the matrix elements of a spherical 
tensor operator can be factorized into two parts. The 
first is a geometric part identified by a Clebsch-Gordon 
coefficient, while the second is a reduced matrix element 
that does not depend on the magnetic quantum numbers. 

To derive the spherical form of any EOM-CCSD 
method, a suitable spherical many-body basis has to be 
defined. The Slater determinant basis where we use the 
particle-hole formalism is substituted with the standard 
form. That means that for a general operator O, we will 
write 

{ab\6\ij) EE (<i>^^|6|$o), (35) 

where the single-particle states labelled a and b are occu- 
pied in the outgoing state, while the single-particle states 
labelled i and j are occupied in the incoming state. All 
single-particle states shared between the incoming and 
outgoing many-body state are dropped from the nota- 
tion. 

In general form, we write a component of the spherical 
basis as 

\a;JaMa), (36) 

where a denotes a particular many-body state, while 
Jc,{Ma) is the total angular momentum(projection) of 
this state. Using the spherical notation, the matrix ele- 
ments of the excitation operator arc written 

r^(J„, Jp) = (a; JMHiW; JpM^), (37) 



where we have dropped the cumbersome sub- and super- 
scripts on the excitation operator in favour of standard 
tensor notation. The matrix elements of the matrix vec- 
tor product in Eq. (|22l) are written 

(a; J„Af,| (RRii) ^ JpMfi) = 

io{a;Jc.Ma\Rij\P;JpMp). (38) 

Now, the Wigner-Eckart theorem allows a factorization 
of the matrix elements into two factors 

CMi/;M>-^J-\\{M')jP;Jp) = 

ujciXMMJc.m'mjp). (39) 

Here, C^^^^j^j^ is a Clebsch-Gordon coefficient and the 
double bars denote reduced matrix elements and does not 
depend on any of the projection quantum numbers. This 
equation is simplified by dividing by the Clebsch-Gordon 
coefficient. This means that for each set of a, /3, Jq, and 
Jp, where J, Jq and Jp satisfy the triangular condition, 
there are (2 J-l-1) x (2 J^-l-l) x (2 J^-|-l) identical equations 
for a given J. Of course, we only need one to solve the 
eigenvalue problem, resulting in a significantly smaller 
dimensionality. We are left with an eigenvalue problem 
where the unknown components of the eigenvectors are 
the reduced matrix elements of the excitation operator 

(a; JqII (hR'')^||/3; J^) ^ uo{a-Jo,\\k^\\[3-Jfi). (40) 

The eigenvalue problem in Eq. (l40l) is the spherical for- 
mulation of the general EOM-CC diagonalization prob- 
lem. For a given excitation operator, both the connected 
operator product and the reduced amplitudes must be 
defined explicitly. 



III. DERIVATION OF SPHERICAL 
2PA-EOM-CCSD 

In this work we derive the spherical formulation of the 
2PA-E0M-CCSD T| method, where the excitation oper- 
ator in Eq. (PT|) has been truncated at the 3P-1H level. 
We write 

R = ^ E + ^ E (41) 

ah abci 

where we have dispersed with the cumbersome sub- and 
subscripts in the operator. 

The unknown amplitudes r are the matrix elements of 
R and defined by 

^ab ^ ($ah|j^|^^^ ^ (a6|R|0) 

^abc ^ ($f ^|R|$„) ^ (a6c|R|*), (42) 
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while we define a shorthand form of the components of 
the matrix- vector product by 



(M) -(<i>'^''|(HR)^|$o) 

(hr)^ =(<i>^ 



HR) |$o> 

c 



(43) 
(44) 



The eigenvalue problem in Eq. ^F2\ can now be posed 
in terms of the individual equations 



(M) 



'HR 



ah 



(45) 



In the spherical formulation, the excitation operator is 
a spherical tensor operator of rank J and projection M , 



ah 

^ ^ ^ {.Ji Jabci Jat 



ahci 

'Jab'^abc 



(8) hiji) 



M 



(46) 



Here the a\{ja) and ai{ji) are spherical tensor operators 
of rank ja and ji respectively, where the latter is the time- 
reversed operator of a\{ji). Standard tensor notation has 
been used to define the spherical tensor couplings. The 
reduced amplitudes are now the reduced matrix elements 
of the spherical excitation operator (|46l) . We write 



(47) 



where ja and jh are coupled to J in left to right order. 
Also 

r";^" (J, Jabc, Jab) = {abc;jajb; Jabjc] J abcW'R^'^ W', ji) , 

(48) 

where we have coupled the angular momentums ja and 
jh to Jab, as well as Jab and jc to Jafcc, also in left to 
right order. The shorthand form of the reduced matrix 
elements of the connected operator product is defined 
analogously by 

(hR'^)"' (J) = {ab-jajb-, J\\ (M'^)^ ||0) (49) 



and 



/ _ ^ \ o,bc 

(^HR j (J, Jabc, Jab) 



{ahc;jajb] Jabjc] JabcW (hR"'^ ^ IN; Ji)- (50) 



The transformations that connect the reduced matrix el- 
ements of R'^ with the uncoupled matrix elements are 
given in Eqs. (|Al6|) - (|Al9j) . 

The spherical eigenvalue problem (j40p can now be writ- 
ten 

(nR'y\j) = Lur''\j) 

(M^) (J, Jabc, Jab) = ^rf ^( J, Jabc, Jab), (51) 



where the amplitudes are now the reduced matrix ele- 
ments defined above. 

In Table |T] we list the main result of this section. The 
first column lists all possible diagrams that contribute to 
the matrix- vector product in Eqs. (|45l) . and (|5T|) . The 
other two columns contain the closed form expressions 
for these diagrams, in both the uncoupled and in the 
spherical representation. All matrix elements and ampli- 
tudes are defined in appendix |Al while the permutation 
operators P(a, 6) and P(a&, c) are defined in appendix IB] 
Note that in the spherical representation the permutation 
operators also change the coupling order. 

The last two diagrams contain the three-body parts 
of the similarity-transformed Hamiltonian (|S]) and have 
been joined in the spherical representation. The details 
of this merger and the definition of the intermediate op- 
erator X are contained in appendix ICl 

Let us briefly discuss how we arrive at the spherical 
diagram expressions, where we will use the first diagram 
in Table H] as an example. This diagram contributes to 
the 2P matrix elements defined in Eq. (|43p in the uncou- 
pled representation and to the reduced matrix elements 
defined in Eq. (|49l) in the spherical representation. 

The first step is to use the transformation in Eq. (jAl6P 
to write the reduced matrix elements (j49l) in terms of the 
uncoupled matrix elements (|43)). This gives us 



(My (J) 



1 

J2 



E 

AlnialTij, 



C: 



'jajbJ 



M 



HR 



ah 



(52) 



where we introduce the notation j — \J2J + 1. The dia- 
gram contributions to the uncoupled matrix elements are 
given by 



(hr)"'' ^ P(a&)H^r'' 



(53) 



where we let the arrow indicate that this is only one of 
several contributions to this matrix element. Here, Hg is 
a matrix element of the one-body part of the similarity- 
transformed Hamiltonian ([8]) and both H^ and r*^"^ are in 
the uncoupled representation. 

Second, we insert Eq. (f53| into Eq. ([52]) and get 



HR 



{J) 



(54) 



Note that for the moment, we have ignored the permu- 
tation operator P(a6) that is a part of the diagram. We 
will discuss this later. 
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Diagram 



Uncoupled expression 



Coupled expression 



P{ab)m{j,)r''%J) 



iP(a6)H 



bm aef 
ef ' m 



xR:y{J,f)r'^\j,f,J,ft„J) 




P(a,6c)H^fr'' 



P{ab,c)mrt 





iP(a6,c)H:>f^ 



1 fjatc e/ 
2 J^e/i' 



P{ab,c)m{j,)rf''{Jab,Jat„J) 



^^P{ab,c)nff{Jab)rf''{Jab.. Jaba, J) 



fJab je Jabe 
jc Jmc jm 
Jabc ji J 



(_l)l+ic+i™-JjJ^ 



jc jm J 1 c 



ji Jabc Jal 



Xm, {J)^im(Jab) 



Table I: All diagrams for the 2PA-EOM-CCSD method with both ordinary and reduced amplitudes and matrix elements. The 
reduced amplitudes and matrix elements are defined in appendix [A] while Xi{J) is defined in appendix [C] Note that the two 
last diagrams are combined into one expression in the spherical formulation and that repeated indices are summed over. 



Third, we use the reverse transformations in 
Eqs. (jAllI) and (IA17P to transform the uncoupled ma- 
trix elements of and r"*^ to the corresponding reduced 
matrix elements. This gives us 

M 

^ 5Z '^Jt Je'^™6.™eC^fj^jj^A,/C;^^^„^j(,f , (55) 

where S is the Kronecker delta and results from the ap- 
plication of the Wigner-Eckart theorem to the matrix 
element of H. Orthogonality of the Clebsch-Gordon co- 



efficients gives 

marrib 

and we are left with 

(HR)"'(J)^H^U)r-(J). (57) 

Note that X]j\/l = 2J+l and that repeated indices are 
summed over. 

Initially, we left out the permutation operator P(a, b) 
that is needed to generate antisymmetric amplitudes. In 
the uncoupled representation this operator is defined as 

P(a&) = i - Pa,b, (58) 
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where 1 is the identity operator and Pa,b changes the 
order of the two indices a and b, but leaves the cou- 
pling order unchanged. Let us apply this operator to 

(My (J). The resuh is 

P(a,6) (llRy\j) ^ 

(nRy\j)-{ba;j,,n;J\\ (hr)^||0), (59) 

where the last matrix element has the wrong coupling 
order compared to the reduced amplitudes defined in 
Eq. (gHl) where 

(M)'" (J) = {ba;jbja;J\\ (M)^ ||0). (60) 

To change the coupling order, we exploit the symmetry 
properties of the Clebsch-Gordon coefficients and write 

{ba;jajb;J\\ (Sr)^I|0> - 

{-iy^+''-'{ba;jbJa;J\\ (SR)^I|0) 

= (_l)Ja+J->-./(^HR)''"(J). (61) 

As a consequence of this, we define the permutation 
operator in the spherical representation to also change 
the coupling order. This results in the following definition 

P{ab) = i - i-iy a,t- (62) 

The total contribution from the first diagram in Table|T] 
in the spherical representation is 

(HR)"'(J)^P(a6)H^(jb)r-(J), (63) 

where P{ab) is defined by Eq. 

We define the more complicated three-body permu- 
tation operators P(a6, c) in the same manner, but they 
must change the coupling order of three angular momen- 
tums. The details have been left to appendix IB] 



IV. RESULTS 
A. Model space and interaction 

All calculations in this section have been done in a 
spherical Hartree-Fock basis, based on harmonic oscilla- 
tor single-particle wave functions. These are identified 
with the set of quantum numbers {nlj} for both protons 
and neutrons, where n represents the number of nodes, 
I represents the orbital momentum, and finally j is the 
total angular momentum of the single-particle wave func- 
tion. 



We will identify the size of the model space by the 
variable 

7Vn,ax = niax(Ar), (64) 

where N = 2n + l, so that the number of harmonic oscil- 
lator shells is iVniax + 1. All single-particle states with 

271 + 1 < TVmax 

are included and no additional restrictions are made on 
the allowed configurations. Thus, iVmax completely de- 
termines the computational size and complexity of the 
calculations. 



max 


Size 


Elements 


Memory 


10 


132 


145 623 788 


1.1 Gb 


12 


182 


587 531302 


4.4 Gb 


14 


240 


1 963 734 704 


14.6 Gb 


16 


306 


5 687 352 954 


42.4 Gb 


18 


380 


14 715 230 212 


109 Gb 


20 


458 


33 622 665 364 


250 Gb 



Table II: The first column contains the size of the single- 
particle basis employed for different model spaces labelled by 
-^max(See text for details). Column two and three list the 
number of matrix elements for the different model spaces and 
the memory footprint of the interaction in our implementa- 
tion. All numbers are based on the coupled representation, 
also known as the jj-scheme basis. 

In Table |ll] we list the size of the single-particle space 
for different values of A'max in the spherical representa- 
tion. In addition, the table includes the total number of 
matrix elements defining the interaction in Eq. ([1]), as 
well as the memory footprint in our calculations. Given 
the memory requirements alone, it is clear that a dis- 
tributed storage scheme is needed from the beginning. 

In addition to the interaction elements, we also need 
to store the Arnoldi vectors in the diagonalization pro- 
cedure. Typically 150 iterations are performed, where 
one vector has to be stored for each iteration. Table IIIII 
lists the size of a single vector for selected target states 
in various model spaces. For a double precision calcu- 
lation, where each element requires 8 bytes of storage, 
the Arnoldi diagonalization would require « 76 Gb of 
memory for a calculation of the J'^ = S"*" state of ^Li 
with A^max = 16. Thus the Arnoldi procedure quickly 
becomes the largest memory consumer in this method. 
In general, there is a steep performance penalty in in- 
creasing the total angular momentum of the target state, 
comparable to increasing the size of the model space. 

Our interaction is derived from chiral perturbation the- 
ory at next-to-next-to-next-to-leading order(N'^LO) us- 
ing interaction matrix elements of Entem and Machleidt 
21j. The matrix elements of this interaction employs a 
cutoff A = 500 MeV and all partial waves up to relative 
angular momentum J^ei — 6 are included. We do not in- 
clude relevant three- and four-body interactions defined 
by the chiral expansion at this order. 
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State 


iVmax = 10 


iVmax = 12 


Afmax = 14 


A^max = 16 


iVmax = 18 


A^max = 20 


''He(0+) 


516 048 


1 323 972 


2 981 930 


6 088 376 


11 513 088 


20 176 104 


^He(l-) 


1 507 930 


3 894 028 


8 808 688 


18 040 354 


34 190 482 


60 Oil 982 


'^He(2+) 


2 391 692 


6 251 128 


14 255 896 


29 364 090 


55 885 624 


98 356 664 


«Li(0+) 


775 992 


1 989 508 


4 478 936 


9 142 216 


17 284 308 


30 285 212 


''Li(l+) 


2 268 746 


5 853 534 


13 234 004 


27 093 632 


51 335 514 


90 080 136 


'^Li(2+) 


3 595 384 


9 391 650 


21 409 878 


44 088 456 


83 893 672 


147 629 532 


«Li(3+) 


4 676 372 


12 438 258 


28 699 916 


59 604 726 


114 125 048 


201 657 602 


i»O(0+) 


1 908 474 


5 022 710 


11 485 808 


23 680 034 


45 071 990 


79 331 610 


i*o(i-) 


5 594 899 


14 802 528 


33 974 801 


70 231 288 


133 940 727 


236 049 974 




8 891 923 


23 794 936 


55 036 119 


114 391 274 


219 038 683 


387 077 788 


1*0(2-) 


8 897 760 


23 803 219 


55 047 530 


114 406 595 


219 058 796 


387 083 193 


1*0(3+) 


11 613 562 


31 596 862 


73 906 056 


154 840 950 


298 237 942 


529 098 382 


1*0(3-) 


11 621 868 


31 608 838 


73 922 708 


154 863 424 


298 267 524 


529 107 862 


1*0(4+) 


13 629 562 


37 905 214 


89 982 332 


190 504 054 


369 757 342 


659 327 780 


i*F(0+) 


2 868 568 


7 545 420 


17 248 686 


35 552 756 


67 658 660 


119 071 548 


^F(l+) 


8 403 602 


22 228 738 


51 009 366 


105 427 688 


201 040 066 


354 285 892 


i*F(2+) 


13 362 878 


35 742 012 


82 642 970 


171 734 254 


328 788 766 


580 957 010 


i*F(3+) 


17 451 568 


47 458 334 


110 973 350 


232 452 890 


447 659 068 


794 095 862 


i*F(4+) 


20 479 376 


56 930 198 


135 106 850 


285 982 274 


554 996 372 


989 530 134 


i*F(5+) 


22 363 324 


63 896 228 


154 444 460 


331 158 558 


648 765 300 


1 163 943 530 



Table III: Size of the many-body space in the diagonalization procedure in the Arnoldi algoritlim for all states calculated in 
this work. All numbers are based on the coupled representation, also known as the jj-scheme basis. 

When we discuss the center-of-mass contamination of 
the final solution, we will also need a softer interaction 
where the short range parts are removed via the similarity 
renormalization group transformation (SRG) [1^. We 
will use a cutoff A = 2.0fm-^ for this purpose. 



B. Center of mass contamination 



Recently, Hagen et al. [1^ showed that the coupled- 
cluster wave function separates into an intrinsic part and 
a Gaussian for the center-of-mass coordinate. This means 
that we do not have to worry about the center-of-mass 
contamination of physical states, as the results will con- 
verge to a separated solution. 

In the EOM-CC approach we make further approxi- 
mations by truncating the many-body basis before the 
diagonalization is performed. It is therefore not clear 
that the final wave functions separates in the same way 
as the coupled-cluster reference state. 

Following the procedure in Hagen et al. [1^ , we calcu- 
lated the expectation value of the center-of-mass energy 
at the proper frequency, for the ground state of ^He, the 
first excited J'^ = 1" state in ^He and the first excited 
J'^ 3+ state in *^Li. The calculations were performed 
in a model space defined by A^max = 14, with an SRG 
interaction with a momentum cutoff at A = 2.0 fm"^. 

Fig. [T] shows the calculated center-of-mass frequency 
as a function of the oscillator parameter Huj defining the 
single-particle basis, while Fig. [5] shows the center-of- 
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^0 15 20 25 30 



Figure 1: Center-of-mass frequency as calculated by the de- 
scription in Hagen et al. [2^ as a function of the oscillator pa- 
rameter fkjj. Three different states are shown - the — 0+ 
ground state of ®He, the = 1^ excited state in ''He and 
the = 3+ excited state in ''Li. 



mass energy at this new frequency. 

Immediately, two things are evident. First, the J"^ = 
f- state in ^He does not separate. The center-of-mass 
energy is very high and we conclude that this state is 
unphysical. There are no documented — 1^ low- 
lying states in the ®He spectrum. Second, the ground 
state of ^He separates, as the center-of-mass energy is 
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— «He(0+) 

— 6He(i-) 

— 6Li(3+) " 









15 



20 



25 



30 



Figure 2: The expectation value of the center-of-mass energy 
is calculated at the center-of-mass frequency to for the — 
0+ ground state of ^He, the = 1~ excited state in ®He and 
the — 3^ excited state in ^Li. 



approximately zero. At low frequencies, the calculation 
is not quite converged, yielding a higher center-of-mass 
admixture in the wave function. 

Let us comment further on the — S"*" state in ^Li. 
Here, the center-of-mass energy is not quite zero, al- 
though very close at approximately 300 keV. This is a 
physical state, as ^Li has this state low in its spectrum. It 
must be said that the Hellman-Feynmann theorem used 
to calculate these values is not exact, so small errors can 
occur. In addition, the calculation might not be fully 
converged, but we saw no significant changes between 
calculations done at A^max = 12 and Amax = 14. Finally, 
it can be that the result is not converged with respect to 
the truncations made in the theory. If this is the case, 
it might be possible to use the center-of-mass contami- 
nation as a measure of convergence. We will touch upon 
this subject in the next two sections, but a detailed anal- 
ysis is left to future work. 



C. Applications to "Li and "He 

For any given reference nucleus, there are only three 
nuclei accessible to the 2PA-E0AI-CC method. Using 
■*IIe as the reference, we can add two protons to calcu- 
late properties of ®Be, two neutrons for ^He and finally 
a proton and a neutron to calculate properties of ^Li. Of 
these, only ^Li and ^He are stable with respect to nu- 
cleon emission and we will only concentrate our efforts 
on these. Looking at the structure of ^Li and ^He, they 
are very different. This is an important point, because 
the quality of the current level of approximation will in- 
evitably depend on the structure of the nucleus. 

^Li is well bound and has four bound states below 
the nucleon emission threshold at 4.433 MeV [l^]. The 
ground state has spin parity assignment J'^ = 1+, while 



the first excited states have J'^ = 3+, 2+ and 0+. The 
J'^ = O"*" ground state in ^He is a two neutron halo state, 
bound by only 800 keV [l^] compared to *He. There are 
no bound excited states - only a narrow resonance at 
1.710 MeV is documented (H- 

First we discuss convergence with respect to the size 
of the model space. Our approach is based on a finite 
basis expansion, where A^max effectively determines the 
numerical cutoff. We increase the cutoff until the correc- 
tions are so small that the uncertainties in our calculation 
outweighs these differences. 



^max = 


10 


^max = 


12 


^max = 


14 


^max = 


16 




Figure 3: Ground state energy of Li as a function of the oscil- 
lator parameter to and the size of the model space A'^max (see 
text for details). 

Fig.[3]shows the calculated total binding energy of ®Li, 
as a function of the oscillator frequency huj. Different 
lines correspond to different model spaces. At A^max = 
16, there is a shallow minimum around huj = 24 MeV and 
in a 10 MeV range including this minimum, the binding 
energy varies by approximately 100 keV. This is less than 
half a percent of the total energy. At low frequencies, the 
energy deviates substantially from the minimum, due to 
the lack of resolution in the single-particle space. 

Note that the gain in binding energy when going from 
A^max = 14 to A^max = 16 is also Very small, about 40 keV. 
We conclude that the binding energy of ^Li is converged 
with respect to the size of the model space (Amax)- The 
energy at fia; = 24 MeV will be our tabulated result. 

The picture is largely identical for the binding en- 
ergy of ^He, only the minimum in energy occurs at 
hoj — 20 MeV. Here the difference in energy between 
the two largest model space is about 140 keV. 

To show that the excited states follow the same pattern 
of convergence as the ground state, in Fig. 0] we plot the 
total energy of the — state in ^Li. As before, the 
energy is plotted as a function of huj and different lines 
correspond to different values of A'max- Our estimation 
of excitation energies will also be taken at specific values 
of huj, where we simply take the difference in energy to 
the ground state at the minimum value. For the calcula- 



11 



' — > ^max = 


= 10 


■ — ■ ^max = 


= 12 


• — • ^max = 


= 14 


"-^ ^max = 


= 16 




Figure 4: The total energy of the first excited — S"*" state 
of ®Li as a function of the oscillator parameter f'lio and the 
size of the model space A'max (see text for details). 
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Figure 6: Excitation energy for the first excited = 2"*" state 
of of ®He, as a function of the size of the model space defined 
by Nmax- The rightmost column shows the experimental value 
from Tilley et al. [H]. 



tion of excitation energies in this section, this causes no 
concern as the minima in energies for all excited states 
are located at the same fko value as the minima for the 
ground states. (See next section for a discussion of cases 
where this is not correct.) 




Figure 5: Excitation energy for selected states of ®Li, as a 
function of the size of the model space defined by A'max. The 
rightmost column shows the experimental values from Tilley 
et al. [13]. 

In Fig. [5] we show the convergence of the excitation en- 
ergies for selected states in the spectrum of ^Li. The hori- 
zontal axis denotes the size of the model space, where the 
rightmost column are the experimental values [l^. All 
excitation energies have been calculated at fiw = 24 MeV, 
which correspond to the minimum of the ground state en- 
ergy. All three states show very little dependence on the 
size of the model space at iVmax = 16. None of the states 
show any significant center-of-mass contamination, so we 
conclude that they are converged. We also found an ad- 
ditional J'^ = 1"*" state higher in the spectrum, but this 



state had a large center-of-mass admixture in the wave 
function and is therefore not included. 

In Fig. [6] we show an equivalent plot for the first — 
2+ excited state in ^He, which is also converged with 
respect to the size of the model space. No significant 
center-of-mass contamination was found neither in this 
state nor the ground state, but this was not the case 
for the J'^ = 1^ state, as already discussed. Note that 
the scale of the vertical axis of this plot is significantly 
smaller than that of Fig.[S]and that all excitation energies 
for ^He were calculated at hoj — 20 MeV. 

Wc also look at some properties of the wave function. 
Although it is not an observable, other expectation values 
might be more sensitive to changes in the wave function 
than the energy. For this purpose we define two proper- 
ties. First, the partial norms are defined by 



ii2p0h)^Y.(^J+l){r^'iJ)y 



(65) 



ab 



n(3pl/i)= i'2Jabc + l){rf'iJab,Jabc,J))\ (66) 

abci 

JabJabc 

where n{2p0h)+n{3plh) = 1. The amplitudes r°'^{J) and 
rf^'^iJab, Jabc, J) are the spherical amplitudes defined in 
Eqs. (|A16P and (|A18p . respectively, while Jx are angu- 
lar momentum labels. Note that the angular momentum 
factors are included so that the partial norms are consis- 
tent between the coupled and uncoupled scheme. These 
norms quantify the part of the wave function in 2P-0H 
and 3P-1H configurations, respectively. 
Second, we define the total weights 



ab 
pw 



(67) 



where the label pw identifies the partial wave content of 
the weight. The sum is over all configurations with this 
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partial wave content, because the weights of individual 
configurations are not stable. Only at high frequencies do 
single configurations dominate. The total weights how- 
ever, are reasonably stable already at this level. 



State 


n{2p0h) 


Dominant configuration(s) 


Weight (s) 


»Li(l + ) 


0.77 


(7rp3/2l^Pl/2), {P3/2y 


0.40, 0.30 


''Li(3+) 


0.76 




0.67 


''Li(0+) 


0.72 




0.60, 0.11 


«Li(2+) 


0.77 




0.640 


'=Hc(0+) 


0.71 


{■Pzlif, (Pi/2)^ 


0.59, 0.11 


'^Hc(2+) 


0.74 


(P3/2)^, (P3/2P1/2) 


0.47, 0.25 


«He(l-) 


0.65 


(S1/2P3/2), (^5/2^3/2) 


0.32, 0.18 





Expt. 


N^LO 
(A = SOOMeV) 


NCSM[25] 


Egs(l ) 


1 nno 
— oi.99o 


— 28.44(0 ) 


— zo.o(o ) 


Eg.(l+) 
Ex (3+) 
Ex(0+) 
Ex (2+) 


0.0 
+2.186 

+3.562 
+4.312 


0.0 

+3.60(4) 
+3.76(3) 
+4.44(4) 


0.0 

+2.91(3) 
+3.30(10) 
+4.10(15) 


«He 


Egs(0+) 


-29.270 


-25.51(14) 


-26.2(5) 


Egs(0+) 
Ex(2+) 


0.0 
1.797 


0.0 

2.35(4) 


0.0 



Table IV: This table shows the 2P-0H partial norms (|65|) . as 
well as the dominant configurations for calculated states in 
both ^Li and ''He. The weights are calculated according to 
Eq. H67p where all nodes for a given partial wave contribute 
to the sum. The numbers in parenthesis are the differences 
of the maximum and minimum values for &j = 20 — 40 and 
indicate the level of convergence of the wave function. 

We list the partial norms and dominant weights of se- 
lected states in ^Li and ^He in Table IIVI The ground 
state of ^Li is dominated by configurations that have a 
proton in the orbital and a neutron in the or 
the orbitals. Note that both partial norms and the 
total weights are more sensitive to changes in the model 
space than the energy and that the weights are becoming 
more unstable approaching the nucleon emission thresh- 
old. Again the J'^ = 1^ stands out with a low 2P-0H 
content. The 2P-0H structure of all the physical states 
in this table conform with the shell-model picture of two 
nucleons in the p-shell, but all states have an additional 
3P-1H content that cannot be neglected. Additional re- 
finements of our method using also 4P-2H configurations 
are necessary to analyze this further. 

In Table |V] we show results with estimated numerical 
uncertainties for the ground and selected excited states 
of both ^He and ®Li. For comparison, both experimental 
values and results from a no-core shell-model (NCSM) 
calculation (2^ are tabulated where data is available. 
Note that the results from the NCSM calculation are 
based on the same interaction as the results from this 
work, but the interaction is renormalized using the pro- 
cedure defined in Suzuki and Lee [1^, before the diago- 
nalization was performed. In addition, the final results 
were extrapolated to an infinite model space. 

Let us discuss the uncertainties indicated by the paren- 
thesis in the table. For the results from this work, listed 
in the second column, the numbers in parenthesis sim- 
ple give the difference in energy between the two lar gest 
model spaces. The results from Navratil and Caurier [25 
give the extrapolation errors in the last column, while the 
experimental energies [24] in the first column are listed 
without uncertainties. 

Fig. [7] shows a graphical representation of the data 
in Table El Compared to the results from the NCSM 
calculation, our results are quite promising. First, the 
ground state energy of ^Li is well within the uncertain- 



Table V: Binding energies and excitation energies for selected 
states of ®Li and ®He with estimated numerical uncertain- 
ties. The uncertainties in our results are the differences be- 
tween the values for the two largest model spaces. The data 
is compared to NCSM results [23| in the rightmost column, 
where the parenthesis list extrapolation errors and experi- 
mental data in the leftmost column. All experimental data 
are from Tilley et al. [24<] . 




^ <^ 



/ 



Figure 7: Excitation levels of selected states in Li and 
calculated using 2PA-EOM-CC(this work) and NCSM 
while compared to experimental values. 



'He, 



ties of the "exact" result, while the ground state energy 
of ^He is just outside. The difference between the two 
nuclei can be explained by the extended spatial distri- 
bution of ^He. Additional correlations are necessary to 
account for this structure. Although the alpha core in 
^He is expected to stay largely unchanged when adding 
two neutrons, the distribution of these extra neutrons are 
biased in one direction. This results in a skewed center- 
of-mass compared to the center-of-mass of the alpha core 
alone. Additional correlations are necessary to absorb 
the resulting oscillations of the alpha core with respect 
to the combined center-of-mass. The spatial distribution 
of ^Li is tighter, so this effect is not that prominent in 
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this nucleus. Second, the ordering of excited states in ^Li 
is reproduced. 

We do however, consistently overestimate the excita- 
tion energies. While the difference for the first J'^ = 0+ 
and 2+ states are small enough to be attributed to the 
difference between the calculations, this is not so for the 
J'^ = 3+ state. Our results differ by about 700 keV for 
this state. Neither the partial norms nor the total weights 
listed in Table ITVl provide any hint of explanation for this 
discrepancy. About 76 percent of the wave function is in 
2P-0H configurations, which is comparable to the ground 
state. The wave function is dominated by configurations 
where both nucleons are in orbitals, which is con- 
sistent with the shell-model picture. Furthermore, the 
level of convergence for this state is no different from the 
other excited states. We do however, find a slight center- 
of-mass admixture in the wave function of the J'^ = 3+ 
state, which is not present in the other states. This was 
illustrated in Fig. [2J where we used a softer renormal- 
ized interaction. Now, the calculated center-of-mass con- 
tamination is too small to be conclusive, especially since 
the calculation is only approximate and the analysis has 
been performed using a softer renormalized interaction. 
We include it because it might indicate that additional 
correlations are needed in the calculation, either in the 
reference or the EOM operator. This matter needs to be 
investigated further, but currently the implementation 
will not allow model spaces large enough for a converged 
description of the center-of-mass admixture in the final 
wave function. 

Using the in-medium similarity renormalization 
group (IM-SRG), Tsukiyama et al. \2t\ performed a sim- 
ilar study with a softer interaction. Here, the = 3+ 
state in ^Li is reproduced on the same level of accuracy as 
the other bound states. We therefore conclude that the 
discrepancy in our result is due to the approximations 
done in the 2PA-E0M-CCSD method. 

Let us also comment on the differences between our 
results and the experimental values. First, we see that all 
excitation energies are overestimated compared to data. 
Again, the J'^ = 3+ states is exceptional, but this has 
been discussed in detail by Navratil and Caurier (25| . 
The matter was resolved by the inclusion of three-nucleon 
forces [28j . which also brought the binding energy very 
close to data. 

For the = 2+ resonance in ^He, this seems not to 
be the case. We also calculated this state using a chiral 
interaction with a different cutoflt of 600 MeV and found 
the excitation energy of this state to be largely cutoff in- 
dependent. This was not so for the excited states in ^Li, 
where especially the J"^ = 3"*" state turned out to be very 
cutoff dependent. We can better explain the discrepancy 
in energy of the resonance in ^He, by missing continuum 
degrees of freedom. The current single-particle basis can- 
not handle the description of both bound, resonance and 
continuum states necessary in this case. While we are 
working on the inclusion of the Gamow shell-model [29] in 
the 2PA-E0M-CC method, as has been done for coupled- 



cluster in Hagen et al. [30j a discussion of this is beyond 
the scope of this article. 

Summing up this section, we conclude that for well 
bound states, with simple structure, the current approxi- 
mation will yield total energies comparable to exact diag- 
onalization. For weakly bound states, continuum degrees 
of freedom are necessary for a quantitative comparison 
to exact results. In addition, we find that we are able 
to reach sufficiently large model-spaces for our results to 
be converged for six nucleons, but for certain states, es- 
pecially the J"" = 3+ state in *^Li, the effects of 4P-2H 
configurations need to be investigated. 



D. Applications to ^^O and ^^F 

When ^^O is used as a reference state, the three iso- 
bars 1^0, i^F and ^^Ne are reachable by the 2PA-E0M- 
CCSD method. All have well bound ground states and 
a rich spectra of bound excited states below their re- 
spective nucleon emission thresholds. The spectrum of 
^*F is especially rich, as the exclusion principle does not 
affect the placement of nucleons in the sd-shell. The 
proton-neutron interaction is responsible for the com- 
pressed spectrum in the fluorine isotope, while strong 
pairing effects in ^^O and ^^Ne result in lower ground 
state energies. In both ^^0 and ^®Ne the spectra are 
opened up and the first excited states are higher in en- 
ergy. Our current focus is on convergence and the viabil- 
ity of this method. Thus, ^^Ne is not explicitly discussed, 
as results are similar to those of ^^O. 
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Figure 8: The ground state energy of as a function of the 
oscillator parameter, fku. Different lines correspond to differ- 
ent model spaces, parametrized by the variable A'^max (|64p . 

Let us first look at the convergence of the binding 
energy of ^^O. In Fig. [SI we plot the ground state en- 
ergy of ^^O as a function of the oscillator parameter 
huj. The different lines correspond to different model 
spaces, parametrized by the variable A^^^ax (EH)- A shal- 
low minimum develops around hu = 32 MeV, where the 
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energy is converged with respect to the size of the model 
space. The difference in energy is about 20 keV when 
the size of the model space is increased from iVmax — 14 
to TVmax = 16. We also see that for a wide range of val- 
ues around the minimum, the ground state energy shows 
very little dependence on the huj parameter. This state 
is exceptionally well converged and is a testament to the 
bound nature of this state. 

A similar result is obtained for the ground state energy 
of ^^F, where the difference in energy is around 160 keV 
between the two largest model spaces. This is almost 
an order of magnitude larger than for the ground state 
of ^^0, but is still well within one percent of the total 
energy. 
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Figure 9: The total energy of the — 3+ state in as 
a function of the oscillator parameter, fku. Different lines 
correspond to different model spaces, parametrized by the 
variable A''max (l64l). 



In Fig. ini we plot the total energy of the first excited 
J'^ = state in ^^O for different model spaces. Here we 
see a shallow minimum develop at huj — 28 MeV. Also 
this state is very well converged, with a difference in en- 
ergy of only about 25 keV between calculations in the two 
largest model spaces. We clearly see that the rate of con- 
vergence is different for different values of hu. This will 
be an issue when we want to calculate excitation ener- 
gies, as different choices of Huj lead to different results. To 
evaluate the excitation energy, we highlight two options. 
We can treat the total energies as variational results and 
choose the lowest energy for the ground state and the 
lowest energy for the J'^ ~ 3+ excited state. Thus, at 
A'max — 16 we calculate the excitation energy as 



£;o+(32MeV), 



(68) 



where i?3+(28MeV) is the total energy of the J'^ — 
3+ excited state, calculated at fiw = 28 MeV, while 
i?Q+ (32MeV), is the ground state energy calculated at 
hJ= 32 MeV. 

The second option is to use the same value of hoj 
for both energies, where we choose the value where the 



ground state has a minimum. Thus, for the current case 
we calculate the excitation energy as 

E^{3+) = £;3+(32MeV) - £:o+(32MeV). (69) 

The difference in energy between these two options is 
minimal if sufficiently large model spaces are used, but it 
will have a significant impact on the rate of convergence. 
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Figure 10: The excitation energies of the first = 2 , 3 
and 4+ excited states in where the best values for huj 
has been chosen for each state. The different columns rep- 
resent different model spaces parametrized by the variable 
A''max (|64p . The rightmost column contains experimental val- 
ues from Tilley et al. [3j|. 
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Figure 11: The excitation energies of the first =2"*", 3"*" 
and 4+ excited states in ^**0 for &j = 32 MeV. The different 
columns represent different model spaces parametrized by the 
variable A^max (|64|) . The rightmost column contains experi- 
mental values from Tilley et al. [3l|. 

The effect is best viewed in Figs. [TOl and ITl] where we 
plot the excitation energy of the first = 2+, S"*" and 4+ 
excited states in ^^0 as functions of the size of the model 
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space. In Fig. [TOl the excitation energy is calculated 
according to Eq. while it is calculated according to 
Eq. (|69|) in Fig. [TTJ For ^^O, this choice will only affect 
the — "i^ state, as the other states all have minimum 
values at huj — 32 MeV. The effect is significant, but we 
feel the first approach of Fig. [10] correctly depicts the level 
of convergence of the excited states. All three states in 
this plot are very well converged, but the J'^ = 4+ state 
displays a slightly larger model space dependence than 
the other two states. 

To further evaluate these states, we look at some of 
the properties of the wave function. First, we checked 
that the wave function separates approximately into an 
intrinsic part and a center-of-mass part, following the 
recipe outlined in Hagen et al. '23'] . None of these states 
showed any center-of-mass contamination. Second, we 
looked at the partial norm defined in Eq. (|65l) and also 
the relative weights of the different configurations in the 
final wave function. In Table IVTl we list these properties 
for all states calculated in this section. 

The J'" = 2+ and 4+ states are dominated by configu- 
rations where two neutrons are in the /2 orbital, but as 
can be seen in Table EH the J"" = 2+ state shows signifi- 
cant admixture of configurations where one neutron is in 
the orbital while the other is in the Si/2 orbital. The 
J'^ = 3+ state is dominated by the latter configuration. 
Thus, all three states conform with the shell-model pic- 
ture of two neutrons in the sd-shell outside an ^^O core. 
They also show the same level of 2P-0H content - around 
70 percent. 

In contrast, we can look at the three negative parity 
states that show up in the spectrum. First, from 
Table I VII we see that the 2P-0II content of the wave 
functions are significantly lower for the negative parity 
states. 
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Figure 12: The excitation energies of the first J'^ = 1 ,2 
and 3" excited states in ^**0 for fiu = 32 MeV. The differ- 
ent columns represents different model spaces parametrized 
by the variable Vmax (|64p . The rightmost column contains 
experimental values from Tilley et al. [3lj . 



Second, Fig. [12] shows that the excitation energies of 
these states are not converged at N^^-^ = 16, and far 
from experimental data. Finally, the center-of-mass con- 
tamination of all three states is rather large and similar 
to the J'^ = 1^ state in ^He (see Fig. [2]). We must there- 
fore conclude that the current truncation is not sufficient 
to reproduce these states. It will be necessary to include 
4P-2II correlations to get the negative parity states right. 
We can understand this by examining how negative par- 
ity states can be produced in ^^0. First, it can be pro- 
duced by placing one neutron in the sd-shell, while the 
other is placed in the p/-shell. If this was the dominant 
configuration, the current truncation would have been 
enough. Secondly, it can also be produced by placing 
two neutrons in the sc?-shell and excite nucleons from the 
p-shell up to the sd-shell. If this excitation is comparable 
in energy to the first configuration, 3P-1H configurations 
start to dominate and 4P-2II configurations are probably 
necessary. 

As a side note, we also see that the wave function of 
these negative parity states are spread over a large num- 
ber of configurations with significant contributions from 
higher partial waves. 

In the spectrum of ^^O, there are three bound J'^ = 0+ 
and 2+ states. The second J'^ = 0+ state is es pec ially 
interesting for this method, as it is a 4P-2H state [S^]. In 
the shell-model language, it is an intruder state, because 
configurations outside the sd-shell are important to get 
this state right. As the current implementation only in- 
cludes 3P-1II configuration, this state can provide clues 
as to what behaviour we can expect from states that are 
not converged with respect to the level of approximation. 
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Figure 13: The excitation energies of the first few = 0+ 
and 2-\- states in for hw — 32 MeV. The different columns 
represents different model spaces parametrized by the vari- 
able A^max (|64p . The rightmost column contains experimental 
values from Tilley et al. [sj]. 

In Table [VI] we list three J"" = 0+ states in ^^O and 
none of them stands out. They all have similar partial 
norms of around 70 percent and are dominated by two 



16 



State 


n{2p0h) 


Dominant configuration(s) 


Weight (s) 




''*O(0+) 


0.69 


(d5/2)^ (Sl/2)' 


0.52,0.12 






0.69 


(Sl/2)" 


0.62 






0.70 


\ o/ / 


0.62 




"0(2^) 


0.70 


(c?5/2)^, (Si/24/2) 


0.37, 0.27 




i«0(2+) 


0.71 




0.44, 0.26 




«0(2+) 


0.70 


(^5/2)^3/2, (51/2^3/2) 


0.49, 0.17 




"0(3^') 


0.71 


(51/2^5/2) 


0.71 






0.71 




0.64 




i»o(ir) 


0.62 


(P3/2rf5/2),(P3/2Sl/2), (^5/2/7/2) 


0.22, 0.17, 


0.12 


"0(2^) 


0.65 


(^3/2^5/2) 


0.49 






0.62 


(^3/2^5/2), (^1/20(5/2), (Sl/2/7/2), (^5/2/7/2) 


0.19, 0.12, 


0.11, 0.10 




0.69 


(^5/2)' 


0.51 




i»F(l+) 


0.73 




0.28, 0.22, 


0.16 


i»F(2+) 


0.72 


{Sl/2d5/2),{d5/2d3/2), (51/2^3/2) 


0.34, 0.20, 


0.14 


i»F(3+) 


0.71 


(^1/2^5/2), (4/2)^ 


0.47, 0.18 




isf(4+) 


0.72 


(4/2)' 


0.64 




i«F(5+) 


0.72 


(rf5/2)' 


0.70 





Table VI: This table shows the 2P-0H partial nornis H65p . as well as the dominant configurations for calculated states in both 
"O and "F. The weights are calculated according to Eq. (|67|) where all nodes for a given partial wave contribute to the sum. 



neutrons in S1/2 and ^3/2 respectively. In Fig. [T3l 

we show the convergence patterns for these states, along 
with those of the three J'^ = 2+ states. All states show 
similar level of convergence and it is not possible to single 
out one of the states. However, if we look at the center- 
of-mass contamination of these states, the third J'^ — 
0+ shows massive contamination, while the other two 
show almost none. We therefore chose to associate our 
J'^ = 0+ state at 11.6 MeV, with the experimentally 
observed J"" = 0+ state at 3.63 MeV. Our ,r = 0+ state 
at 3.4 MeV is closer in energy, but including effects of 
three- nucleon forces, pushes this state higher in energy 
and very close to the experimentally observed = 0+ 
state at 5.34 MeV. A similar effect occurs among the 
J'^ = 2+ excited states, but it is less prominent. Here, 
the center-of-mass contamination were negligible for all 
three states. 

Let us summarize our discussion of missing many-body 
correlations. We have identified four different markers to 
indicate missing physics. Unfortunately, none of them 
can be used quantitatively and all must be evaluated si- 
multaneously to form a general picture. First, the par- 
tial norm can be used to differentiate between different 
states. From these calculations it seems that a 2P-0H 
norm of around 70 percent is the standard. A lower par- 
tial norm, might indicate the need for 4P-2H correlations. 
Second, we look at the total weights in the different con- 
figurations. If the weights are spread over many differ- 
ent configurations, as they were in the wave functions of 
the negative parity states calculated in ^^O, this might 
also indicate missing 4P-2H configurations. Third, we 
look at the convergence patterns and if energies converge 



slowly, this probably means that something is missing 
from the calculation. In weakly bound states for exam- 
ple, continuum effects result in the need for additional 
resolution in the single-particle basis. Finally, we look at 
the level of center-of-mass contamination present in the 
wave function. None of these arguments can be analyzed 
in detail, before 4P-2H configurations are included. We 
are working on this, but computationally, it will only be 
possible to include these configurations in a small single- 
particle space. If the 3P-1H and 4P-2H configurations 
are defined only in a so-called active space around the 
Fermi level, the computational cost might be manage- 
able. This has been done successfully in Gour et al. [3j| 
and should prove to be a valuable approximation also in 
this method. The formation of a correlated alpha cluster 
around the Fermi level is important in this mass region 
and can hopefully be accounted for using a minimal set 
of 4P-2H configurations. 

Let us also look at the convergence of selected states 
in ^*F. In Fig. [T3] we show the excitation energy of the 
first few states in ^*F for different model spaces. Here all 
states are relative well converged, with only the J'^ = 4+ 
state showing some model space dependence at A^max = 
16. None of these states have significant center-of-mass 
contamination and all partial norms are on the same level 
as can be seen in Table IVII This table also shows that 
the dominant weights are spread over several different 
configurations, but all within the sc?-shell. According to 
our earlier discussion, we can find no fault with these 
states. 

In Fig. [TSl we plot selected states from the excitation 
spectra of ^^O, ^^F and ^*Ne. Only states that are con- 
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Figure 14; The excitation energies of the first — , 
2^, 3^, 4^ and 5^ states in '^^F, calculated at fto; = 
30 MeV. The different columns represents different model 
spaces parametrized by the variable A'max (|64|l . The right- 
most column contains experimental values from Tilley et al. 
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-0- n- 



18q 
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N^LO 
(A = 500MeV) 


Egs(0+) 


-139.807 


-130.00(2) 


Egs(Or) 
Ex(2+) 
Ex (4+) 
Ex(Oj) 
Ex (2+) 
Ex (2+) 
Ex(0+) 
Ex (3+) 


0.0 
+1.982 
+3.554 
+3.633 
+3.92 
+5.255 
+5.336 
+5.3778 


0.0 

+1.50(1) 
+2.72(7) 

NA 
+3.38(2) 
+7.27(18) 
+3.40(15) 
+4.12(1) 


18p 


Egs(l+) 


-137.370 


-129.75(16) 


Eg.(l+) 
Ex (3+) 
Ex(0+) 
Ex(5+) 
Ex (2+) 
Ex (4+) 


0.0 
+0.937 
+1.041 
+1.121 
+2.523 
+4.652 


0.0 

+2.15(4) 
+2.92(4) 
+2.77(2) 
+3.43(1) 
+5.81(12) 


i^Ne 


Egs(0+) 


-132.143 


-122.56(15) 


Egs(0+) 
Ex (2+) 
Ex (4+) 
Ex(3+) 


0.0 
+1.887 
+3.376 
+4.561 


0.0 

+1.43(3) 
+2.59(6) 
+3.83(4) 



Table VII: Ground and excitated state energies for the A = 18 
system, compared to experimental data. The experimental 
energies are from Tilley et al. [3ll]. The number in parenthe- 
sis indicates the level of convergence and is the difference in 
energy between calculations in the two largest model spaces. 



is much more open compared to experiment than that of 
^^0. The level ordering is largely reproduced, except for 
the cases already discussed. A more detailed comparison 
to experimental data in this region, is difficult to perform 
before both three-nucleon forces and continuum degrees 
of freedom are included [1, |38| - |40| . 



Figure 15: Excitation energies for selected states in O, F 
and ^®Ne compared to experimental values from Tilley et al. 




sidered good are plotted and compared to their experi- 
mental values. 

Table I VIII lists the numerical values, used in Fig. [12] 
together with the ground state energies for future com- 
parison. The difference in energy between calculations in 
the two largest model spaces, is listed as the uncertainty 
of our results in this table. The experimental values are 
from Tilley et al. (sj 

It is clear that the spectrum in ^^O and ^*Ne is com- 
pressed compared to data, while for ^®F, the spectra is 
more open. This is consistent with other calculations us- 
ing different interactions [33 - [37j . where the ^^F spectra 



V. CONCLUSIONS AND OUTLOOK 

We have presented the spherical version of the 2PA- 
EOM-CCSD method, appropriate for the calculation of 
energy eigenstates in nuclei that can be described as two 
particles attached to a closed (sub-)shell reference. The 
method has been evaluated in both A = 6 and A = 
18 nuclei, where the results were mostly converged with 
respect to the single-particle basis. 

In comparison with a full diagonalization, both ground 
state and excited state energies were in general very ac- 
curate. However, a singular excited state in ^Li deviated 
significantly from the "exact" result, showing the need to 
include additional correlations like 4P-2H configurations 
for the accurate treatment of complex states. For sim- 
ple states, where a 2P structure is dominant, the current 
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level of truncation is adequate. 

In comparison with experiment, both three-nucleon 



forces and a correct treatment of the scattering contin- 
uum are needed to refine our results. 
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Appendix A: Reduced matrix elements 



The reduced matrix elements of a spherical tensor op- 
erator T^^ of rank J and projection M, are defined ac- 



cording to the Wigner-Eckart theorem, 



{a-Jc.M^\Tij\p-JpMp) 



^MMaMa. 



(a; J„||T^ 



(Al) 

Here, a and /3 are general labels representing all quan- 
tum numbers except angular momentum and its projec- 
tion, while Ja{Ma) and JplAlfj) are the total angular mo- 
mentum(projection) of the bra and ket state respectively. 
The double bars denote reduced matrix elements and 
does not depend on any of the angular momentum pro- 
jections and C'ljlj^^j^j is a Clebsch- Gordon coefhcient. 

In coupled-cluster, the unknown amplitudes are the 
matrix elements of the cluster operator Tq 



it = {a;jama\'TQ\i;jimi) 
tfj = {ab; jamajbmb\To\ij; jzmijjTTij) 



where the operator sub- and superscript identifies the 
cluster operator as a scalar under rotation. Also, the 
labels abij . . ., denote single-particle states and we have 
singled out the angular momentum(projection) in the la- 
bels jairria) etc. 

Now, the reduced matrix elements of the cluster oper- 
ator are defined according to Eq. (|A1[) . For example. 



{ah; jamajbmb\T^\ij:,jimijjmj 



C: 



r 



a J a J a 



JabMa 

JijMi. 



manihMab m-irrijMij OMijMab 



X {ab;jajb; Jab\\'i°\\ij;jijj;Ji3), 



where the sum and the first two Clebsch-Gordon coeffi- 
cients come from the coupling of ja and jt to Jab and 
of ji and jj to Jtj, in that specific order This expres- 
sion is simplified by the explicit evaluation of the third 
Clebsch-Gordon coefficient 

where 5a, p is the Kronecker delta. We get 
(a5; jamajbrnb\Tl\ij-,jirnijjmj) = 



EfljajbJ 
^niarribM 

JM 



C^'^^M(«^;jaJ;,;J||TO|Kj;j.j,;J). 



When this specific coupling order is used(left to right) 
and no confusion will arise, we will use a shorthand no- 
tation for the reduced matrix elements, defined by 

t^(J) = (a;J||T°||z;J) (A2) 
fi^iJ) = {ab-jajb-, ^||T°||*j; J). (A3) 



The transformations between the reduced amplitudes 
and the original amplitudes of the cluster operator are 
defined as 



t^ = S. 



/-ijajbJ f-i3i3jJ >ab 



rriimj M 



t: 



ab 



M''ij h 



(A4) 
(A5) 

(A6) 
(A7) 



JM 



where we use the convention J — \/2J + 1. 

As the similarity-transformed Hamiltonian(H) is a 
scalar under rotation as well, the shorthand form of its 
reduced matrix elements are defined analogously 

H^(J) = (p;J||H"||g;J) (A8) 
KKJ) = {pq;jpjg;Jm"\\rs;jrjs;J). (A9) 

The transformations between the original and reduced 
matrix elements are given by 



H?(J) = <5,„,,=jH^ 



(AlO) 
(All) 



1 



HP? 



JM 



(A13) 



The original matrix elements of the similarity- 
transformed Hamiltonian are defined in Jansen et al. [7]. 

The r-amplitudes are the matrix elements of the exci- 
tation operator in Eq. pS)) 



„ab 



„abc 



{ab; jamajbmb\Rii\0) 

{abc; jamajbmbjcmclRijli; jinii), 



where the operator is now a general tensor operator of 
rank J and projection M. The bra side contains up to 
three indices and we have to couple three angular momen- 
tum vectors to be able to define the reduced amplitudes. 
Using rf'"^ as an example, we couple from left to right 
and get 



{abc; jamajbmbjcmc\RM\i; jiirii 



E 

JabMcb 
JabcMabc 



C: 



jajbJab 



c 



JabjcJaba 



c 



malTLijAIab ^'^ab'l^c ^^abc ^'^"l^i^^ab, 



{abc;jajb; Jabjc, Jabc\\'R.'^\\i;ji 



where the last Clebsch-Gordon coefficient is due to the 
Wigner-Eckart theorem. We let the order of the angular 
momentum labels on the bra-side specify the coupling 
order, where ja and jb couples to Jab- In turn. Jab and jc 
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couples to Jabc- When this coupling order has been used 
and no confusion will arise, we will use the shorthand 
notation for the reduced elements, defined by 

r'^^J) = {ab-]aJb;jW\\o) (A14) 

rf" {J, Jabc, Jab) = {abc;jajb\Jabjc\Jabc\\^''\\i]ji)- 

(A15) 

In the shorthand notation, the transformations be- 
tween the reduced and the original amplitudes is 



ah _ i^dajbJ ab( j\ 

abc( T T ^ T — ^ \^ abcf^jajbJab 

'i {'J, 'Jabc, 'Jab) — ^ 'i '-"m^rribA' 



(A16) 

(A17) 



"^abc MMab<:Mab 
rriaTribmcmi 

„ f-^JabjcJaba f^JjiJabc 



(A18) 



„ahc 



J2 rf^iJ^ Jabc, Jab) 

Jabc^abc 

JabM^b 

V riiaibMab r^JabOaJaba /~<JjiJabc 

^ ^rriarribMab Mabm^Maba ^MrmMa. 



(A19) 



Appendix B: Permutation operators 



The diagrams in TableUcontain permutation operators 
that guarantees an antisymmetric final wave function. In 
the uncoupled formalism, these were simple and defined 

by 



F{ab) = i - Pa,fc 
P(ab,c) = i-Pa.,c-Pb, 



(Bl) 
(B2) 



where Pa,b permutes indices a and b and 1 is the identity 
operator. In the spherical formalism, this simple form is 
not adequate, as a specific coupling order is used in all 
reduced amplitudes. To see why, let us apply Pa.b to a 
reduced amplitude r'^''{J) 

Pa.b{^b;jajb;Jm'm = {ba;jajb;Jm'm 

While this coupling order is the correct order when cal- 
culating the contribution to (^S-E^) , the reduced am- 
plitudes are defined in a different coupling order. To 
compensate, we change the coupling order and introduce 
a phase(See Edmunds [20*1 for details) 

\ba;jaJb;J) = i-iy-+''-'\ba;jbja;J). 

Thus, we define 

P(a6) = i-(-iy»+^^-^Pa,6, (B3) 



and applied to the reduced amplitude r°*{J) this has the 
correct form 

P{ahy\j) = r''\j) - (-1)J"+J'''-V°(J). 

The permutation operators Pa,c and Pb,c in Eq. (jB2p 
are a bit more complicated, as they involve three particle 
states. We use standard expressions for coupling three 
angular momentum 



\cba\ jajb; Jabjc, JabcMabc) 



^ ^ '^cb'Jab 



cb 



jc jb J, 
ja Jabc Jab 



X \cba] jcjb] J cbj a] JabcMabc) 
\acb- jaib-, J abic, JabcMabc) = ^1 '^"^''^^^ 

jc ja Ja 



\jb Jabc Jab 
X \acb] jajc, Jacjb] JabcMabc 



Thus, we define 



P{ab, c) = i -f ^ JcbJab 



ja Jabc Jab 



Pa.. 



'y^^(^_iyb+ja~.Jab + 'Jac 



jc jb Jcb 
c '-^a 

jc ja J a 



'Jab'Jac ^ 1 ^' T J { ^b.c, 

Jb 'Jabc 'Jab 



(B4) 



Since this anti symmetrization contributes a significant 
part of the overall calculation, all diagrams containing 
this operator is applied only once to the sum of all dia- 
grams containing this operator. 



Appendix C: Three-body parts of H 

There are two three-body matrix elements of H that 
contributes to the 3P-1II amplitudes. But since the orig- 
inal Hamiltonian does not contain three-body elements, 
these deserve special attention. These elements can be 
factorizcd in just the same way as the coupled-cluster 
amplitude equations, to reduce the computational cost 
of these diagrams. The two three-body contributions to 
the 3P-1H amplitudes are 



1 

1 



efi m 



1 

^P(a,6c)Hr/e,C/. 
These terms are factorized to get 



--P(a, bc)ll^Jr^ft'^, + -P{a, bc)R':pr^h'„Tr. 

nn 

2' 



1 



(HR)f 

where we have defined the intermediate 



-F(a,6c)Cx^ 



ef ''n 



(CI) 
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Note that we have swapped the indices to facihtate the 
angular momentum couphng. Now the angular momen- 
tum coupling in the coupled-cluster amplitudes match 
the coupling in the 3P-1H amplitudes, so we do not need 
to break these couplings when rewriting the diagram in 
a spherical basis. 

In the spherical formulation, it is clear that Xm ^^e 
reduced matrix elements of the tensor operator x' , which 
has the same rank as R'-'. This is a consequence of the 
scalar character of H. By coupling the matrix elements 
in Eq. (jCip to form reduced matrix elements, we get the 



following expression for the reduced matrix elements of 

+ E (-i)-+-+^--^^|f f\ 

Jett/c -^efc -^'^^ -^^n 

y<rlf^{J,f,J,f,,J)W^p{J,}). (C2) 



